'''
Created on Dec 9, 2010

@author: oabalbin
'''
import exome.gatk_cluster.cluster_jobs_header as jh

def snps_call_SomaticSnipper(ref_genome, normal_bam_file, tumor_bam_file, snps_calls_file, path_to_snipper):
    '''
    '''
    
    snipper_command=path_to_snipper+'bam-somaticsniper'
    args = [snipper_command, '-q',str(1), 'processSomatic']+tumor_bam_file+\
           normal_bam_file+[ snps_calls_file]

    args= [a.replace(',',';') for a in args]
    command = ",".join(args).replace(',',' ').replace(';',',')

    return command


def main_snniper_caller(cohort_dict, run_dict):
    '''
    '''
    ref_genome=run_dict['ref_genome']
    path_to_snipper=run_dict['path_to_snipper']
    short_wtime=run_dict['short_wtime']
    extra_mem=run_dict['use_mem']
    out_dir = run_dict['out_dir']
    node_processors = run_dict['num_cores']
    single_processor=1
    my_email = run_dict['myemail']
    
    jobn='sss'
    i=0
    for cohort, cohort_samples in cohort_dict.iteritems():
        # Note that in this case there list_bam_files_normal
        # correspond to a merged bam file of normals or tumors
        # for that cohort
        jobn=jobn+str(i)
        list_bam_files_normal = cohort_samples['normals']
        list_bam_files_tumor = cohort_samples['tumors']
        
        snps_calls_file = out_dir+cohort+'.snipper_snps.'
        command = snps_call_SomaticSnipper(ref_genome, list_bam_files_normal, list_bam_files_tumor, snps_calls_file, path_to_snipper)
        
        jobidvc = jh.qsub(jobn+'_vc', command, single_processor, cwd=None, walltime=short_wtime, pmem=extra_mem, 
                           deps=None, stdout=None, email_addresses=my_email)
                
        # Look how to convert the final varscan output into a vcf format
        i+=1
        return jobidvc
    
        
if __name__ == '__main__':
             
    cohort_dict = {'Prost':{'normals':['/nobackup/med-mctp/oabalbin/test/Prost.normals.sam_snps.mpileup.vcf'],
                            'tumors':['/nobackup/med-mctp/oabalbin/test/Prost.tumors.sam_snps.mpileup.vcf']}}


    run_dict = {'path_to_gatk':'/nobackup/med-mctp/sw/bioinfo/gatk/GenomeAnalysisTK-1.0.4705/',
                 'resources_folder':'/nobackup/med-mctp/sw/bioinfo/gatk/GenomeAnalysisTK-1.0.4705/resources/', 
                 'rscipt_path':'/home/software/rhel5/R/2.10.1-gcc/bin/Rscript',
                 'path_to_picard':'/nobackup/med-mctp/sw/bioinfo/picard/picard-tools-1.35/', 
                 'path_to_bwa':'/nobackup/med-mctp/sw/bioinfo/alignment/bwa/bwa-0.5.8a/',
                 'path_to_sam':'/nobackup/med-mctp/sw/bioinfo/samtools/samtools-0.1.10/',
                 'path_to_vcftools':'/nobackup/med-mctp/sw/bioinfo/vcftools/',
                 'path_to_varscan':'/nobackup/med-mctp/sw/bioinfo/varscan',
                 'path_to_snipper':'/nobackup/med-mctp/sw/bioinfo/SomaticSnipper/',
                 'use_mem':8000, 'num_cores':6,
                 'nmismatches':2,
                 'short_wtime':'24:00:00',
                 'long_wtime':'60:00:00',
                 'ref_genome':'/nobackup/med-mctp/sw/alignment_indexes/gatk/hg19/hg19.fa',
                 'ref_genome_bwa':'/nobackup/med-mctp/sw/alignment_indexes/bwa/hg19/hg19.fa',
                 'snpdb_file':'/nobackup/med-mctp/sw/alignment_indexes/gatk/hg19/dbsnp132_00-All_processed.vcf',
                 'indeldb_file':'/nobackup/med-mctp/sw/alignment_indexes/gatk/hg19/dbsnp132_00-All_processed.vcf',
                 'target_exons':'/nobackup/med-mctp/sw/alignment_indexes/agilent/hg19/exome/Agilent_SureSelect_All_Exon_G3362.v2.hg19.bed',
                 'recal_analysis_outputdir':'/nobackup/med-mctp/oabalbin/test/recal_analysis/',
                 'temp_dir':'/nobackup/med-mctp/oabalbin/test/temp/',
                 'qsubfile':'/nobackup/med-mctp/oabalbin/test/',
                 'out_dir':'/nobackup/med-mctp/oabalbin/test/',
                 'myemail':['alebalbin@gmail.com']
                 }
    
    main_snniper_caller(cohort_dict, run_dict)
